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Solutions of the time-dependent, mass-averaged Navier-Stokes equations are compared in detail with ex- 
perimental results obtained on an axisymmetric “bump” model at a transonic Mach number that produced an 
extensive separated flow region. In addition, an inverse boundary method is evaluated for this type of flow. The 
Cebeci-Smith algebraic and the Wilcox-Rubesin two-equation turbulence models used in the Navier-Stokes 
calculations both predict the maximum boundary-layer displacement thickness generated by the interaction 
reasonably well, with the details of the flow best described with the two-equation formulation. However, both 
models predict a shock location substantially farther aft on the bump than observed experimentally. This error 
in shock location was slightly less with the two-equation model (0.12 chord compared with 0.16 chord). In the 
vicinity of the shock, the calculations predict a more rapid increase in turbulent shear stress than observed in the 
experimental results; this more rapid increase is believed to be the cause of the poor predictions in shock 
position. 


Nomenclature 

c ~ chord of bump 

Cj = skin-friction coefficient, T w fVip^u 2 m 

k = kinetic energy of turbulence, (u' 2 + t/ 2 + w' ; )/2 
L = length-scale function in two-equation turbulence 
model 

M = Mach number 

p - static pressure 

p t = total pressure 

Re = Reynolds number 

R f = turbulence Reynolds number 
t = thickness of bump 

u = velocity in x direction 

v = velocity in y direction 

x = axial coordinate parallel to model axis 

y - distance from model surface in direction per- 

pendicular to model axis 
- no rmali zed distance from model surface, 
yyfr^p w /n w 

&(x) = nondimensional pressure gradient parameter 

6 = boundary-layer thickness 

6* = displacement thickness 

7 * = exponential damping term in two-equation tur- 

bulence model 

p = molecular viscosity 

H t = turbulent eddy viscosity 

p = fluid density 

r , = turbulent shear stress 

t w = wall shear stress 

to = turbulent dissipation rate 

Subscripts 

e = boundary-layer edge conditions 

w = wall conditions 

oo = freestream conditions 
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Superscripts 

()' = fluctuating quantity 

( A = time-averaged quantity 

( “ ) = mass-averaged quantity 

Introduction 

I N the past several years, considerable advances have been 
made in prediction of pressure distributions and, 
accordingly, lift performance for lifting surfaces at transonic 
conditions. However, accuracy of these methods deteriorates 
as the interaction between the shock wave and turbulent 
boundary layer on the suction side of the lifting surface 
strengthens with increasing angle of attack or Mach number. 
In general, reasonably accurate predictions can be obtained if 
the interaction is weak enough that boundary-layer separation 
does not occur. When separation does occur, large 
displacement effects can be introduced by the viscous flow to 
produce a pronounced change in the resultant pressure 
distribution. Hence, these viscous displacement effects must 
be accurately predicted if a reasonably accurate surface 
pressure distribution is to be obtained. Moreover, viscous 
flow behavior must be accurately described if the correct 
shock position is to be obtained, since the shock position is 
closely coupled to the viscous flow under conditions of 
separation. 

The inaccuracies in predicted surface pressure distributions 
for separated cases have been generally attributed to 
inadequate modeling of the turbulent Reynolds stresses. M 
There also is some concern about whether the computational 
grids used in these solutions were sufficiently fine to 
adequately resolve the shock wave/boundary-layer in- 
teraction regions. 

A lack of comprehensive experimental data has made it 
difficult to draw any firm conclusions about either the nature 
of the inadequacies in the computational methods or the way 
the predictions could be improved. Some recent experiments 
have been conducted on airfoils in which turbulent Reynolds 
stress data were obtained using laser velocimetry . 5 7 These 
data, however, were not obtained in the immediate vicinity of 
the separation point in any of the experiments because the 
boundary layer was extremely thin at separation. 
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Recently, an experiment 8 was conducted on an axisym- 
metric-flow model designed to simulate the flow about an 
airfoil at transonic conditions. This model was developed 
specifically to study in detail the flow development in the 
immediate vicinity of flow separation induced by a shock 
wave. The model consisted of an annular circular-arc bump 
affixed to a thin-walled cylinder aligned with the oncoming 
flow. The cylinder extended well ahead of the bump to allow 
development of a relatively thick turbulent boundary layer at 
the leading edge of the bump. The axisymmetric con- 
figuration provided a flow that was free of three-dimensional 
effects and facilitated acquisition of near-surface measure- 
ments by the laser velocimeter technique. Detailed 
measurements of the mean velocities (streamwise and vertical) 
and the turbulent Reynolds stresses from upstream of 
separation through reattachment were obtained for 
=0.875. 

The objective of this paper is to provide additional insight 
into the prediction of transonic separated flows by making 
detailed comparisons of numerical solutions with the experi- 
mental results of Ref. 8. Two numerical methods are used in 
these comparisons: 1) solution of the mass-averaged Navier- 
Stokes equations by the method of Ref. 9, and 2) solution of 
the compressible boundary-layer equations by the inverse 
method described in Ref. 10. Because of their general ap- 
plicability to separated flows, the mass-averaged Navier- 
Stokes equations are used primarily in these comparisons. The 
inverse boundary-layer solutions are included because of their 
potential advantages over Navier-Stokes equation solutions in 
turbulent flow analysis; their most obvious advantage is a 
substantial saving in computational time. 

Solutions to the mass-averaged Navier-Stokes equations are 
presented for two turbulence models: the algebraic Cebeci- 
Smith 11 eddy-viscosity model and the Wilcox-Rubesin 12 two- 
equation model. For the inverse boundary-layer method, 
solutions are presented for the Cebeci-Smith eddy-viscosity 
model only. 

Description of the Experiment 

A photograph of the flow model used in the experiment of 
Ref. 8 is shown in Fig. 1 ; flow would be left to right. The thin- 
walled cylinder, 15.2-cm o.d., extended 61 cm upstream of the 
bump leading edge. The bump had a thickness t of 1 .9 cm and 
a chord length of 20.3 cm (2 t/c= 19%). Its leading edge was 
joined to the cylinder by a smooth circular arc of radius 18.3 
cm that was tangent to the cylinder at 3.33 cm upstream and 
to the bump at 2.05 cm downstream of the intersection of the 
arc of the bump with the cylinder. The streamwise distance x 
in this paper is referenced to the leading edge of the bump, 
excluding the fairing (i.e., the intersecting point of the arc of 
the bump with the cylinder). 



Fig. 1 Axisymmetric “bump” model. 


The data utilized in this paper were obtained in the Ames 
2 x 2-Foot Transonic Wind Tunnel, which is a closed-return, 
variable-density, continuous-running tunnel with 21% open 
porous-slotted upper and lower walls. Test conditions were a 
freestream Mach number of =0.875 and a unit Reynolds 
number Re/ m = 13.6x 10 6 /m. With the long initial section of 
the model and this unit Reynolds number, a fully developed 
turbulent boundary layer was ensured without im- 
plementation of a transition strip. 

Local instantaneous velocity data were obtained with a two- 
component laser velocimeter system (described in more detail 
in Ref. 6) from which the turbulent Reynolds stresses and 
kinetic energies were realized, as well as the local mean 
velocities. 

The data presented in this paper were reduced using the 
two-dimensional weighting factor suggested in Ref. 13 to 
correct for velocity biasing. In the experiment, the velocity 
fluctuations in the cross-stream direction were not measured. 
Thus, an assumption was necessary to obtain the kinetic 
energy of turbulence, k= {u^ 2 + v' 2 + w’~ 2 )/ 2; it was as- 
sumed that w' 2 = {u’ 2 + v' 2 )/2. 

Numerical Methods 

Navier-Stokes Solutions 

The partial differential equations used to describe the mean 
flowfield were the time-dependent, mass-averaged Navier- 
Stokes equations for axisymmetric flow of a compressible 
fluid. When the two-equation turbulence model was used, 
these equations were augmented by two additional partial 
differential equations: one for the turbulence kinetic energy, 

{u' 2 + v' 2 + w' 2 ) / 2, and another ^for the square of the 
dissipation rate, qj 2 . Modeling constants were necessary in 
both of these equations for closure. 

Computations using the algebraic model were carried out 
with and without the Cebeci-Smith pressure gradient 
correction to the Van Driest damping constant. Only minor 
variations were obtained, and the results presented herein 
were obtained without this correction. A value of A + =26 
was used for the Van Driest constant. 

In the case of the Wilcox-Rubesin two-equation model, the 
turbulent eddy viscosity ft, is expressed in terms of k and m: 

M, = Y */>(*/«) (1) 

In Eq. (1), y * is an exponential damping term dependent on a 
turbulent Reynolds number /?, = p ( ^fkL ) //x. This Reynolds 
number is based on a length scale of the turbulence L defined 
as L = yfk/o). Although longitudinal curvature terms are in- 
cluded in the mean equations, no additional curvature terms 
have been added to the model equations. The complete set of 
equations solved, including the equations for the two tur- 
bulence models employed, is described in Ref. 9. 

The numerical procedure used was the basic explicit second- 
order, predictor-corrector, finite-difference, time-splitting 
method of MacCormack, 14 modified by the efficient explicit- 
implicit-characteristic algorithm of Ref. 15. A description of 
the method used, along with its adaptation to multiequation 
turbulence model equations, is contained Refs. 9 and 16. Also 
contained in Refs. 9 and 16 is a description of the boundary 
conditions and of the special procedures used for the tur- 
bulence model variables and equations. 

In the computational domain, a mesh was developed that 
allowed different point spacing in each coordinate direction. 
Normal to the flow direction, an exponentially stretched fine 
point spacing was used to resolve that part of the flow where 
viscous effects were important; the predominantly inviscid 
outer flow was described using a uniform coarse point 
spacing. In the flow direction, a variable-point spacing was 
used, with most of the points concentrated along the aft 
portion and immediately downstream of the circular-arc 
bump. The number of points and the minimum spacing used 
normal to the flow varied with the turbulence model. The 
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distance of the first y mesh point from the model wall was 
selected small enough so that the solutions were independent 
of this spacing. For the algebraic model, the first mesh point 
of the model was taken within the viscous sublayer (typically 
within a minimum>' + *>'V r r^p7/Mw< 4 )i for th ^ l ' vo ' 
equation model, the minimum required was smaller by a 
factor of 10 or 20. Typical mesh sizes were 95 points in the 
streamwise direction and 50 points normal to the model (with 
30 to 35 points in the viscous layer). Transition from the 
exponentially stretched fine mesh to coarse mesh occurred 
near the outer edge of the viscous layer in the separated zone. 
In the flow direction the mesh spacing varied from 0.3 cm 
near the shock location to 5.0 cm far upstream and down- 
stream. 

The boundaries of the computational mesh extended in the 
vertical direction from the model surface to the vicinity of the 
wind-tunnel wall and in the flow direction from x= — 140 to 
90 cm (x = 0 at the bump leading edge). The upstream 
boundary conditions were prescribed by uniform freestream 
conditions. The downstream boundary was positioned far 
enough aft of the interaction region that all the gradients in 
the flow direction could be set to zero. This boundary con- 
dition was verified by moving the location of this downstream 
boundary and observing substantially unchanged numerical 
results. The model surface was impermeable, and no-slip 
boundary conditions were applied with a constant wall 
temperature. Additional details concerning the boundary 
conditions for the two-equation model near the model surface 
are discussed in Ref. 9. At the outer boundary (the wind 
tunnel wall), inviscid solid-wall boundary conditions were 
used. Details concerning this procedure are contained in Ref. 
3. Since the experimental tests were performed in a square test 
section that had porous upper and lower walls, and because 
the code calculates axisymmetric flow, several numerical 
solutions were made for various solid-wall radii to ensure that 
the calculated results were not significantly affected by the 
walls. The calculated position of the shock wave and the 
separated zone remained unchanged when the wall radius was 
varied from 1 .5 to 2.5 chords. The solutions presented are for 
a wall radius at 1.78 chords from the model centerline, 
corresponding to a How equivalent to that of the wind tunnel 


Inverse Boundary-Layer Solutions 

The numerical procedures of Ref. 10 were used to solve the 
two-dimensional compressible boundary-layer equations. The 
equations were solved in either the direct or the inverse mode, 
using the generalized Galerkin method. For the direct mode, 
which is applicable to attached flow, the pressure gradient 
parameter 0(x) was specified as in conventional boundary- 
layer marching solutions. Conversely, for the inverse mode, 
which is applicable to attached or separated conditions, the 
wall shear (or, equivalently, the skin-friction coefficient C f ) 
was specified and 0(x) became part of the solution. The 
method was fourth-order accurate in the direction normal to 
the surface and second-order accurate in the streamwise 
direction. To permit the computational domain to grow in the. 
streamwise direction, the Levy-Lees-Dorodnitsyn trans- 
formation was applied. With this transformation, good 
resolution, normal to the surface, can be obtained near the 
leading edge of a model, even though the boundary layer 

becomes much thicker downstream. 

In regions where the flow is attached, solutions were ob- 
tained in a forward marching manner. Within the separation 
zone, an iterative sweeping procedure employing locally 
downwind differencing was used to account for upstream 
influence. The method was not extended to treat axisymmetric 
flows* however, errors associated with the neglect of trans- 
verse ’curvature were sufficiently small for the present flow 
case to allow an evaluation of the applicability of the method. 

Solutions were obtained for 25 and 35 nodal points normal 
to the surface. Near the surface, the nodal points were spaced 
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Fig. 2 Surface pressure distribution. 


logarithmically. In the outer two-thirds of the boundary layer, 
the points were spaced quasi-uniformly to resolve the wake 
flow. At midchord of the bump, the first nodal point was 
located very near that for the Navier-Stokes, Cebeci-Smith 
solution. With the Levy-Lee-Dorodnitsyn transformation, the 
innermost grid point was approximately seven times closer to 
the surface near the leading edge of the model. No essential 
difference was observed between the 25- and 35-nodal-point 
solutions. 

Comparison of Navier-Stokes Solutions 
with Experiment 

The surface pressure distributions predicted from the 
solution of the mass-averaged Navier-Stokes equations for the 
Cebeci-Smith algebraic and the Wilcox-Rubesin two-equation 
turbulence models are compared with the experimental 
surface pressure distribution in Fig. 2. As shown in Fig. 2, the 
predictions of shock location did not agree well with the 
experiment. The Cebeci-Smith model and the Wilcox-Rubesin 
model predict shock locations 0.16 and 0.12 chord lengths 
downstream of the experiment, respectively. An error of this 
magnitude in an airfoil application would most likely result in 
a gross misprediction of the lift and moment. Included in Fig. 

2 are the locations of the boundary-layer separation and 
reattachment points (both predicted and measured). 

The reason for the differences in the experimental and 
predicted shock locations is believed to be attributable to the 
inadequacies of these two turbulence models rather than to 
any deficiencies in the experimental data or the computational 
mesh used in the numerical solutions. To establish the validity 
of the experimental pressure distribution, further testing of 
this model was performed in a larger wind tunnel (the Ames 
6x6 Foot Supersonic Wind Tunnel) where the influence of 
the tunnel walls would with certainty be negligible. For the 
same test conditions, in this larger facility, the shock location 
was within l^o chord of that observed in the Ames 2x2 Foot 
Wind Tunnel. These tests were performed after the 
publication of Ref. 17. 

To confirm that the y mesh spacings were sufficiently fine 
in the calculations, comparisons were made with a boundary- 
layer solution obtained with the boundary-layer method of 
Ref. 10. The specified pressure distribution was that predicted 
from the mass-averaged Navier-Stokes solution, using the 
Cebeci-Smith turbulence model. In the boundary-layer 
solution, the same Cebeci-Smith turbulence model was used. 
Excellent agreement was obtained between the boundary- 
layer and the Navier-Stokes Cebeci-Smith model solutions 
from near the leading edge of the model to just upstream of 
the predicted shock location, where the boundary-layer 
thickness becomes quite small due to the strong favorable 
pressure gradient. The computed profiles just upstream of the 
predicted shock location are compared in Fig. 3. It is at this 
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point— just upstream of the shock— that the requirements on 
the y mesh spacing are most critical. The computed profiles 
are plotted with a logarithmic scale in y to expand the near- 
wall region. As seen in Fig. 3, excellent agreement with the 
boundary-layer solution, which is fourth-order accurate i ny, 
was obtained near the surface. The small differences near the 
outer edge of the boundary layer are due to the increase in 
static pressure of the potential flow with distance away from 
the model. 

A check was made to see if a finer streamwise mesh spacing 
in the vicinity of the shock would produce any significant 
change in shock position. Second solutions for both tur- 
bulence models were obtained with a streamwise grid spacing 
of 0.15 cm (approximately one-sixth the approaching 
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boundary-layer thickness) in the region of x/c from 0.5 to 0.8; 
the shock position remained essentially unchanged. 

For both streamwise grid spacings, the computational mesh 
was sufficiently fine to produce a curved shock. If the 
streamwise mesh spacing is too coarse, a normal shock will be 
predicted, as demonstrated in Ref. 4. 

Just upstream of the experimental shock location, both 
turbulence models predict mean velocity profiles that com- 
pare quite favorably with the experimental data as seen in Fig. 
4a. The predicted profiles in Fig. 4a are essentially identical 
except between >*0.05 and 0.25 cm. In this interval, the two- 
equation model solution is in better agreement with the ex- 
periment. The fullness of the boundary layer evident in Fig. 
4a is a result of the strong favorable pressure gradient im- 
posed on the boundary layer ahead. In Figs. 4b and 4c, the 
calculated and measured turbule nt sh ear stresses and kinetic 
energies are compared. Note that u f v’ rather than pu' v' is 
used in the turbulent shear stress comparison. This con- 
vention is used throughout the paper because it is the former 
quantity that the laser velocimeter actually measures. Ex- 
perimental results obtained at x/c = 0.563 and 0.625 are 
presented in Fig. 4. Between these two stations, the boundary- 
layer properties changed only slightly in the calculations and 
in the experiment, with the exception of the measured tur- 
bulence properties in the outer part of the boundary layer. 
This feature of the measurements will be discussed later. The 
predicted turbulent shear stress and kinetic energy 
distributions agreed well with the experiment in the region 
very close to the model surface, where there exists a 
significant mean velocity gradient. The maximum turbulent 
shear stresses predicted by the Cebeci-Smith and the Wilcox- 
Rubesin models were -3.2xlO~ 3 and -2.5x10 3 , 
respectively, in this region. Unfortunately, the measurement 
station nearest the wall was not sufficiently close to the model 
surface for the maximum values of the actual flow to be 
determined. 

In the outer region of the boundary layer, significant shear 
stresses and kinetic energies, which were not predicted by the 
calculations, were measured. The fluctuations present in this 
region are believed to be primarily due to nonstationary 
compression waves rather than to turbulence fluctuations. A 
part of the total shear stress in this region may be due to the 
slow decay of the turbulent shear stress originally created by 
the pressure rise near the leading edge of the bump; however, 
the increase in shear stress that occurs between the two 
stations in this outer region of the boundary layer appears to 
be due to nonsteady compression waves, which would have 
the effect of producing a negative w, v correlation. The jitter in 
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Fig. 4 Comparison with experiment upstream of shock: a) mean velocity; b) turbulent shear stress; c) turbulent kinetic energy. 
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Fig. 5 Displacement thickness distribution. 


the compression waves is thought to be a result of the large- 
scale fluctuations within the turbulent boundary layer. It does 
not appear that these stresses play an important role in 
separation since the progression toward separation occurs in 
the near-wall region. When the surface pressure starts to rise, 
the shear stress near the wall quickly grows to a level much 
larger than that shown at the outer edge of the boundary 
layers in Fig. 4b. 

To the degree that boundary-layer conditions upstream of 
the shock have been defined experimentally, and to within the 
experimental accuracy of the measurements, the Wilcox- 
Rubesin model appears to have accurately predicted the 
actual boundary-layer properties in this region. Also, 
although the mean velocity profile just upstream of the shock 
was not predicted quite as well with the Cebeci-Smith model, 
differences between the predicted and experimental results 
were small. From the evidence presented, it would appear that 
the poor predictions of the shock location were not due to an 
inadequate definition of the boundary layer in the near-wall 
region, immediately upstream of the shock. 

To illustrate the global differences in the computed and 
experimental viscous flows, a comparison is made (Fig. 5) of 
displacement thickness 6V At the trailing edge, where the 
displacement thickness reaches its maximum values, in both 
the predictions and the experiment, agreement is reasonably 
good. Notice that the predicted displacement thicknesses 
increase at a significantly faster rate just downstream of the 
shock than was observed in the experiment. Not too much 
significance should be placed on any slight differences in the 
predicted displacement thicknesses shown in Fig. 5, because 
there is a degree of uncertainty in choosing the integration 
limits when the boundary-layer thickness is not well defined. 
Also, there is some uncertainty in the experimental 
displacement thickness and shape factors upstream of the 
shock because of the extremely thin boundary layer in this 
region. 

Development of the mean velocity, turbulent shear stress, 
and turbulent kinetic energy profiles from upstream of the 
experimental separation point through reattachment is shown 
in Figs. 6, 7, and 8, respectively. Oil flow visualizations in- 
dicated separation at x/c^ 0.7. In comparing these results, 
recall that the shock locations for the two numerical solutions 
differed considerably from that of the experiment. Also as a 
consequence, the Mach number just upstream of the shock 
was considerably higher in the calculations than in the ex- 
periment (approximately 1.4 as compared with 1.32). Thus, 
the computed boundary layers incurred an overall pressure 
rise substantially higher than that of the experiment. It is 
instructive to note from the present comparisons that, with 
the stronger shock, the two turbulence models do predict a 
maximum boundary-layer thickness and deficit close to that 
observed in the experiment. In contrast, mass-averaged 
Navier-Stokes calculations 18 for a 64A010 airfoil section that 
has a flat-top suction peak ahead of the shock severely un- 
derpredicted the boundary-layer growth at stall with the shock 
farther aft than experiment. 


As evident from Figs. 6 and 7, the two-equation turbulence 
closure model produced a solution that more closely 
represents the experimental data. For instance, the mean 
velocity profile shape in the separated region for the Wilcox- 
Rubesin model is in much better agreement with the ex- 
periment. Also, the slow decay of the Reynolds shear stress 
downstream of reattachment is qualitatively predicted with 
this model and not with the algebraic model. This result 
differs from the solution presented in Ref. 17, which showed 
the Wilcox-Rubesin model having a more rapid decay in shear 
stress downstream of reattachment than the Cebeci-Smith 
solution. This anomaly was a result of a programming error. 
To prevent program instabilities during the time development 
of the flow, the program artificially sets a maximum turbulent 
kinetic energy as a function of the local boundary-layer 
thickness. For the solution presented in Ref. 17, the program 
calculated an incorrect boundary-layer thickness near reat- 
tachment (an order of magnitude too small), resulting in 
suppressed values of turbulent kinetic energy and shear stress 
in this region. 

The major difference between the calculations and the 
experiment is in the rate at which the boundary-layer 
properties change in the immediate vicinity of the shock. As 
evident in Fig. 5, the displacement thickness in this region 
increases much more rapidly in the calculations. This is also 
true of the maximum turbulent shear stress within the 
boundary layer, as graphically illustrated in Fig. 9. The 
Cebeci-Smith model predicts an almost discontinuous rise in 
turbulent shear stress at the shock, whereas the Wilcox- 
Rubesin model solution shows a somewhat more gradual 
increase as should be expected since that model is theoretically 
formulated to account for nonequilibrium effects in the flow. 
As seen in Fig. 9, however, the predicted rate of growth in the 
turbulent shear stress at the shock is significantly greater than 
that observed in the experiment. Also, the predicted shear 
stress gradients dr^dy , in the inner part of the boundary 
layer, tend to be larger in this region. Evidently, the actual 
flow is not able to generate these larger shear stress gradients, 
which allow the predicted boundary layers to negotiate a more 
rapid rise in pressure. As a result, the shock in the experiment 
is positioned well forward of the predicted shock locations. 

Near separation there is a possibility that the normal stress 
gradients are sufficiently large to produce a substantial in- 
crease in the effective pressure gradient. The Wilcox-Rubesin 
formulation makes provision for this by including an additive 
pressure due to the kinetic energy of the turbulence in the two 
momentum equations. Effects of the kinetic energy of tur- 
bulence on the pressure are not taken into account in the 
Cebeci-Smith model, but this apparently made no significant 
difference. 

In summary, that the more sophisticated two-equation 
formulation does not predict substantially different flow 
behavior than that predicted with the simple algebraic eddy- 
viscosity model is somewhat discouraging. Local equilibrium 
turbulence models, such as the Cebeci-Smith model, in 
general tend to overpredict the turbulent boundary layer’s 
ability to negotiate strong adverse pressure gradients since 
they do not account for the slow response of the larger tur- 
bulent eddies. Thus, the prediction of the shock wave too far 
aft with the Cebeci-Smith model was expected. Flowever, 
since the two-equation formulation attempts to account for 
nonequilibrium effects of the turbulence, substantially better 
agreement with the experiment was expected. The com- 
parisons presented in this paper indicate that the Wilcox- 
Rubesin model as presently formulated does not correctly 
model these effects in the region of the shock. 

On the optimistic side, the maximum displacement 
thickness was predicted quite well with no modifications of 
either turbulence model. This is believed to be one of the few 
cases in which that has been accomplished for a flow with 
extensive separation, without resorting to some extreme 
modification (generally with little or no physical basis) of 
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Fig. 6 Mean velocity comparison. 



EXPERIMENT 

NAVIER STOKES SOLUTIONS 

C S MODEL 

W R MODEL 




established turbulence models. Moreover, where the predicted 
mean velocities are in relatively good agreement with the 
experimental results, so are the predicted turbulent shear 
stresses and kinetic energies. 

Evaluation of Inverse Boundary-Layer Method 

As noted in the introduction, the range of applicability of 
the inverse boundary-layer method to separated flows is not 


well established. One way to test the validity of the boundarj 
layer approximations in separated-flow regions is to compar 
inverse boundary layer results directly with solutions to th 
Navier-Slokes equations. In this study, the Navier-Sloke 
solution with the Cebeci-Smith turbulence model discuss© 
previously was used as the standard. To ensure a consisten 
comparison, this same turbulence model was incorporated i, 
the boundary-layer method. An inverse boundary-laye 
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Fig. 9 Streamwise variation in maximum lurbulenl shear stress. 


solution with exactly the same skin-friction distribution as 
that from the Navier-Stokes solution was desired. However, 
due to some unknown stability problem in the inverse 
boundary-layer method, a solution could not be obtained for 
the exact wall-shear distribution predicted from the solution 
of the time-averaged Navier-Stokes equations. In Fig. 10 the 
wall-shear distribution for which a solution (designated as 
solution A) was obtained is compared with that from the 
Navier-Stokes solution. (Solution B will be discussed in the 
following paragraph.) Overall, this solution agrees quite well 
with that obtained from the Navier-Stokes equations. The 
predicted displacement thicknesses and pressure distributions 
are included in Fig. 10. In Fig. 11, a comparison of mean 
velocity profiles at some selected stations is given. 

These results are extremely encouraging; however, it was 
found that the same pressure distribution shown in Fig. 10 
could be obtained with a different specified wall-shear 
distribution. This nonuniqueness of the solution was 
discovered when the wall-shear iteration procedure, described 
in Ref. 10, was used to obtain convergence to the Navier- 
Stokes solution pressure distribution from an arbitrary initial 
guess of the wall shear. When this approach was taken, a 
totally different flowfield was predicted. To verify that this 
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Fig. 10 Skin-friction, pressure, and displacement thickness 
distributions. 
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did not result from surface pressure being an insensitive 
parameter for use as a convergence criterion, a solution was 
obtained using the predicted pressures of the inverse bound- 

soutlon A as the convergence criteria. A very 
slightly negative wall shear in the separation bubble was used 

«. t'n m. 8U TK and convergence criterion for p/p, was 
** ‘° °-® 01 - T*> e results 'his solution (designated 

solution B), which are included in Fig. 10, clearly demonstrate 
a nonuniqueness in solutions. 

Self-similar solutions of the laminar boundary-layer 
equations do exhibit two solutions for the same value of the 
pressure gradient parameter 0: an attached and a separated 
olution. To the authors’ knowledge, however, two 
solutions for the same pressure distribution have not been 
obtained m any application of the inverse boundary-layer 
method to non-self-similar separated flows. 

In Ref. 17, it was postulated that the multiple-solution 
behavior may have been a consequence of the Van Driest 
damping term, however, subsequent solutions have been 
obtained with the Van Driest damping removed which still 
exhibit this nonuniqueness. From an application viewpoint 

rr: ely ’ of Ihe dis P |acem ent thickness 

distributions of Fig. 10 can be compatible with a solution of 

the outer inviscid flow. Thus, with a coupled inviscid-inverse 
boundary-layer approach, the realization of solution B would 
not be possible. 

The present results suggest that the inverse boundary-layer 
method can be a viable approach for analyzing shock-induced 
separated flows. Of course, it will not be applicable for flows 
where the pressure gradients across the viscosity layer play an 
important role in the How development. For the solutions to 
the mass-averaged Navier-Stokes equations presented here 
this was not apparently the case. 

Concluding Remarks 

It has been attempted in this paper to perform very careful 
mass-averaged Navier-Stokes calculations for a transonic 
shock wave/boundary-layer interaction flow that has been 
extensively documented experimentally and is free of un- 
h^if three-dimensional effects and tunnel wall or 
blockage effects. The goal was to obtain solutions with a 
sufficiently refined computational mesh whereby any dif- 

™mh S f V ? be ! Ween the s °l ut ' or is and the experiment 
could be isolated to the turbulence closure models used in the 
calculations. It is believed that this goal was met and that the 
relatively poor predictions of shock location obtained were 
due to inadequacies in the immediate vicinity of the shock of 
the two turbulence closure models used. Compared with the 
experimental results, both models predict too rapid a rise in 
the turbulent shear stress in this region. Overall, the Wilcox- 


Rubesin two-equation model did better at describing the 
response of the boundary layer to the shock wave than did the 
Cebeci-Smith algebraic model. However, the results of the 
present work indicate that there are still major improvements 
needed in turbulence closure modeling before strong transonic 
shock wave/boundary-layer interactions can be reliablv 
treated theoretically. 3 

A secondary goal of this paper was to evaluate the inverse 
boundary method, which because of the shorter com- 
putational times of the method could prove extremely 
valuable in separated-flow parametric studies. For the present 
now problem, it was found that the inverse boundary-layer 
method could, for all practical purposes, provide the same 
results obtained with the mass-averaged Navier-Stokes 
equations. A nonuniqueness in solutions, however was 
observed of which any user of this method should be aware in 
its application. The undesirable second solution, fortunately 
had a displacement thickness distribution that would be in- 
compatible with an outer inviscid How solution. 
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